  library(foreign)
  library(nnet)
  library(effects)
  
  mnl <-read.dta("gelpi_replication_Rcode.dta")
  timetable <- multinom(disapptimetable4 ~ posevents + negevents + posbush + negbush, data = mnl, subset=mnl$bushapp4==0)
  summary(timetable)
  fit.eff <- Effect("posevents" , timetable, xlevels=list(posevents=c(0,1)),fixed.predictors=list(given.values=c(negevents=0,posbush=0,negbush=0)))
  data.frame(fit.eff$model.matrix, fit.eff$prob, fit.eff$se.prob, fit.eff$lower.prob, fit.eff$upper.prob)
